#pragma once

#include <iostream>
#include <fstream>
#include <iomanip>
#include <set>
#include <map>
#include <array>
#include <cassert>
#include <algorithm>
#include <memory>
#include <cstring>

#include "json/json.h"
#include "base64.h"


namespace TETGEN_EXT {
  template <class ARRAY, class PrEqual, class PrLess>
  void Unique(ARRAY& src, PrEqual pre, PrLess pr, size_t i0 = 0)
  {
    std::sort(src.begin(), src.end(), pr);
    typename ARRAY::iterator iter;
    iter = std::unique(src.begin() + i0, src.end(), pre);
    src.erase(iter, src.end());
  }

  template<class ARRAY>
  void Unique(ARRAY& src)
  {
    std::sort(src.begin(), src.end());
    typename ARRAY::iterator iter;
    iter = std::unique(src.begin(), src.end());
    src.erase(iter, src.end());
  }


  template <class T>
  class Triple : public std::array<int, 3>
  {
  public:
    Triple() {}

    Triple(const T* data)
      :std::array<int, 3>(data)
    {}

    Triple(const T& a, const T& b, const T& c)
    {
      (*this)[0] = a;
      (*this)[1] = b;
      (*this)[2] = c;
    }

    Triple<T>&   set(const T& a, const T& b, const T& c)
    {
      (*this)[0] = a;
      (*this)[1] = b;
      (*this)[2] = c;
      return *this;
    }

    Triple<T>& set(const T* da)
    {
      (*this)[0] = da[0];
      (*this)[1] = da[1];
      (*this)[2] = da[2];
      return *this;
    }

    void  reverse()
    {
      std::swap((*this)[0], (*this)[1]);
    }
  };

  template <class Key>
  struct hash
  {
    size_t operator()(const Key& s) const
    {
      return (size_t)(s.hash());
    }

    bool operator()(const Key& s0, const Key& s1)const
    {
      return (s0 == s1);
    }
  };

  template <class ARRAY>
  int BSearch(const ARRAY& src, const typename ARRAY::value_type& val, size_t i0 = 0)
  {
    int mid = 0;
    int begin = i0, end = src.size() - 1;
    while (begin <= end)
    {
      mid = (end + begin) >> 1;
      if (val < src[mid])end = mid - 1;
      else if (src[mid] < val) begin = mid + 1;
      else return mid;
    }
    return -1;
  }

  template <class ARRAY>
  int BSearch1st(const ARRAY& src, const typename ARRAY::value_type& val, size_t i0 = 0)
  {
    int begin = BSearch<ARRAY>(src, val, i0);
    while (begin > 0)
    {
      if (src[begin - 1] < val ||
        val < src[begin - 1])
      {
        break;
      }
      else
        begin--;
    }
    return begin;
  }

  typedef std::vector<int>    DIntArray;
  typedef std::vector<double> DRealArray;
  typedef std::map<int, int>  Int2IntMap;
  typedef std::map<int, double>   IntRealMap;

  typedef std::array<int, 2>  Int2;
  typedef Triple<int>         Int3;
  typedef std::array<int, 4>  Int4;

  typedef std::vector<Int3>    DInt3Array;

  template<class T>
  /**
  * @brief The KeyIndex3 class
  */
  class  KeyIndex3 : public  std::array<T, 3>
  {
  public:
    KeyIndex3(const T& k0 = T(), const T& k1 = T(), const T& k2 = T())
    {
      set(k0, k1, k2);
    }

    KeyIndex3(const T* k)
    {
      set(k[0], k[1], k[2]);
    }

    KeyIndex3(const  std::array<T, 3>& k)
    {
      set(k[0], k[1], k[2]);
    }

    void  set(const T& k0, const T& k1, const T& k2)
    {
      if (k0 < k1)
      {
        if (k1 < k2)
        {
          (*this)[0] = k0;
          (*this)[1] = k1;
          (*this)[2] = k2;
        }
        else
        {
          (*this)[2] = k1;
          if (k0 < k2)
          {
            (*this)[0] = k0;
            (*this)[1] = k2;
          }
          else
          {
            (*this)[0] = k2;
            (*this)[1] = k0;
          }
        }
      }
      else
      {
        if (k2 < k1)
        {
          (*this)[0] = k2;
          (*this)[1] = k1;
          (*this)[2] = k0;
        }
        else
        {
          (*this)[0] = k1;
          if (k2 < k0)
          {
            (*this)[1] = k2;
            (*this)[2] = k0;
          }
          else
          {
            (*this)[1] = k0;
            (*this)[2] = k2;
          }
        }

      }
    }

    bool operator<(const KeyIndex3<T>& src)const
    {
      if ((*this)[0] < src[0]) return true;
      if (src[0] < (*this)[0]) return false;
      if ((*this)[1] < src[1]) return true;
      if (src[1] < (*this)[1]) return false;
      return (*this)[2] < src[2];
    }

    template<typename T1>
    bool  operator==(const KeyIndex3<T1>& src)const
    {
      if ((*this)[0] == src[0])
      {
        if ((*this)[1] == src[1])
        {
          return  (*this)[2] == src[2];
        }
      }
      return false;
    }

    bool operator!=(const KeyIndex3<T>& src)const
    {
      return !(*this == src);
    }
  };


  inline size_t Hash(size_t h1, size_t h2)
  {
    return size_t(((h1 << 16) | (h1 >> 16)) ^ h2);
  }



  inline size_t Hash(const KeyIndex3<int>& val) { return Hash(Hash(val[0], val[1]), val[2]); }

  inline size_t Hash(const Int3& val)
  {
    return Hash(KeyIndex3<int>(val[0], val[1], val[2]));
  }

  struct HashInt3
  {
    size_t operator()(const Int3& val) const { return Hash(val); }
    bool operator()(const Int3& val0, const Int3& val1)const
    {
      KeyIndex3<int> k0(val0[0], val0[1], val0[2]);
      KeyIndex3<int>  k1(val1[0], val1[1], val1[2]);
      return k0 == k1;
    }
  };

  static const  short TET_FACE_VERTEX[4][6] = { { 1,3,2,8,9,5 },{ 0,2,3,6,9,7 },{ 0,3,1,7,8,4 },{ 0,1,2,4,5,6 } };

  template<class T0, class T1>
  class MyPair : public std::pair<T0, T1>
  {
  public:
    MyPair() {}
    MyPair(const T0& _Val1, const T1& _Val2)
      : std::pair<T0, T1>(_Val1, _Val2)
    {}

    bool operator<(const std::pair<T0, T1>& src)const
    {
      return this->first < src.first;
    }

    bool operator==(const std::pair<T0, T1>& src)const
    {
      return this->first == src.first;
    }
  };

  template<typename ARRAY, typename Pr >
  class ArrayHashProxyC
  {
  public:
    typedef MyPair<size_t, int>       HashItem;
    typedef typename ARRAY::value_type  value_type;
    ArrayHashProxyC(const ARRAY& array_, Pr pr = Pr(), size_t i0 = 0)
      :_array(array_), _pr(pr), _i0(i0)
    {
      (this->_clean) = 0;
    }

    void set_dirty(bool isdirty = true) { if (isdirty) (this->_clean) = 0; else (this->_clean) = 1; }
    int  find(const value_type& it, Int2& intv);
    int  whereIm(int whi)const { return (this->_hash_vals)[whi].second; }
    int  find(const value_type& it);
    int  count(const value_type& it);
    int  count_duplicate(DIntArray&);
    void getDuplicate(DIntArray&);
    void getDuplicate(std::function<void(const DIntArray&)>);
    void group(DIntArray&, DIntArray&);
    const ARRAY& data()const { return _array; }

    const typename ARRAY::value_type& operator[](size_t i)const
    {
      return _array[_hash_vals[i].second];
    }
    void  update();
    template<typename T2, typename HashFun2>
    int  find(const T2& it, HashFun2&, Int2& intv);

    template<typename T2, typename HashFun2>
    int  find(const T2& it, HashFun2&);
  protected:
    void  hash();

    size_t hash_(const value_type&)const;
    static    size_t Clean(const ARRAY& myarray, Pr& pr, std::vector<HashItem>& hash_vals, size_t j0, size_t sz);
    size_t    i0()const { return _i0; }
  protected:
    const ARRAY&            _array;
    Pr                      _pr;
    size_t                  _i0;
    std::vector<HashItem>   _hash_vals;
    int                     _clean;
  };


  template<typename ARRAY, typename Pr = hash<typename ARRAY::value_type> >
  class ArrayHashProxy : public ArrayHashProxyC<ARRAY, Pr>
  {
  public:
    typedef MyPair<size_t, int> HashItem;
    ArrayHashProxy(ARRAY& array_, Pr pr = Pr(), size_t i0 = 0)
      :ArrayHashProxyC<ARRAY, Pr>(array_, pr, i0)
    {}

    void unique(ARRAY* deleted = NULL);
    void erase_duplicate(ARRAY* deleted = NULL, int nErase = 0);
  protected:
    void  cleanup(const DIntArray&, int cnt_new, ARRAY* deleted);
  };


  template <typename ARRAY, typename Pr>
  int ArrayHashProxyC<ARRAY, Pr>::find(const value_type& it)
  {
    Int2 intv;
    return this->find(it, intv);
  }

  template <typename ARRAY, typename Pr>
  template<typename T2, typename HashFun2>
  int  ArrayHashProxyC<ARRAY, Pr>::find(const T2& it, HashFun2& fun)
  {
    Int2 intv;
    return this->find<T2, HashFun2>(it, fun, intv);
  }

  template <typename ARRAY, typename Pr>
  int ArrayHashProxyC<ARRAY, Pr>::count_duplicate(DIntArray& nums)
  {
    int cnt_all = 0;
    if (this->_array.empty())
      return cnt_all;
    this->hash();
    nums.clear();
    size_t j, ji = this->i0(), sz = this->_array.size();
    nums.resize(sz, 0);
    while (ji < (sz))
    {
      size_t j0 = ji;
      size_t posi = (this->_hash_vals)[ji].second;
      nums[posi] = 1;
      ji = Clean(this->_array, (this->_pr), (this->_hash_vals), j0, sz);

      nums[posi] = ji - j0;
      for (j = j0; j < ji; j++)
        nums[(this->_hash_vals)[j].second] = nums[posi];
      if (nums[posi] > 1)
        cnt_all++;
    }
    return cnt_all;
  }

  template <typename ARRAY, typename Pr>
  void ArrayHashProxyC<ARRAY, Pr>::getDuplicate(std::function<void(const DIntArray&)> functor)
  {
    if (this->_array.empty())
      return;
    this->hash();
    size_t j, ji = this->i0(), sz = this->_array.size();
    while (ji < sz)
    {
      size_t j0 = ji;
      size_t posi = (this->_hash_vals)[ji].second;
      ji = Clean(this->_array, (this->_pr), (this->_hash_vals), j0, sz);
      int cnt = ji - j0;
      if (cnt > 1)
      {
        DIntArray rts;
        for (j = j0; j < ji; j++) {
          rts.push_back((this->_hash_vals)[j].second);
        }
        functor(rts);
      }
    }
  }

  template <typename ARRAY, typename Pr>
  void ArrayHashProxyC<ARRAY, Pr>::getDuplicate(DIntArray& rts)
  {
    this->getDuplicate([&rts](const DIntArray& vals) {
      for (auto it : vals)
      {
        rts.push_back(it);
      }
    });
  }


  template <typename ARRAY, typename Pr>
  void ArrayHashProxyC<ARRAY, Pr>::group(DIntArray& groups, DIntArray& nums)
  {
    if (this->_array.empty())
      return;
    this->hash();
    nums.clear();
    size_t j, ji = this->i0(), sz = this->_array.size();

    while (ji < sz)
    {
      size_t j0 = ji;
      ji = Clean(this->_array, (this->_pr), (this->_hash_vals), j0, sz);
      for (j = j0; j < ji; j++)
      {
        groups.push_back((this->_hash_vals)[j].second);
      }
      nums.push_back(groups.size());
    }
    (this->_clean) = 1;
  }

  template <typename ARRAY, typename Pr>
  int ArrayHashProxyC<ARRAY, Pr>::count(const value_type& it)
  {
    Int2 intv;
    this->find(it, intv);
    return intv[1] - intv[0];
  }


  template <typename ARRAY, typename Pr>
  template<typename T2, typename HashFun2>
  int  ArrayHashProxyC<ARRAY, Pr>::find(const T2& it, HashFun2& fun, Int2& intv)
  {
    if (this->_array.empty())
      return -1;
    if ((this->_clean) == 0)
      this->update();

    HashItem data(fun(it), 0);

    intv[0] = BSearch1st((this->_hash_vals), data, this->i0());
    if (intv[0] < this->i0())
      return -1;
    size_t sz = this->_array.size();
    while (intv[0] < sz && data.first == (this->_hash_vals)[intv[0]].first)
    {
      if (fun(this->_array[(this->_hash_vals)[intv[0]].second], it))
      {
        break;
      }
      intv[0]++;
    }

    if (intv[0] >= sz || data.first != _hash_vals[intv[0]].first)
      return -1;
    intv[1] = intv[0] + 1;
    while (intv[1] < sz)
    {
      if (data.first != (this->_hash_vals)[intv[1]].first)
        break;

      if (fun(this->_array[(this->_hash_vals)[intv[1]].second], it))
      {
        intv[1]++;
      }
      else
        break;
    }
    if (intv[0] < _hash_vals.size() && fun(_array[_hash_vals[intv[0]].second], it))
      return _hash_vals[intv[0]].second;
    return -1;
  }

  template <typename ARRAY, typename Pr>
  int ArrayHashProxyC<ARRAY, Pr>::find(const value_type& it, Int2& intv)
  {
    return this->find<value_type, Pr>(it, this->_pr, intv);
  }

  template <typename ARRAY, typename Pr>
  size_t ArrayHashProxyC<ARRAY, Pr>::Clean(const ARRAY& myarray, Pr& pr, std::vector<HashItem>& hash_vals, size_t j0, size_t sz)
  {
    size_t j, j1, gd;
    j1 = j0 + 1;
    for (; j1 < sz; j1++)
    {
      if (hash_vals[j1].first != hash_vals[j0].first)
        break;
    }

    for (j = j0 + 1; j < j1; j++)
    {
      if (!pr(myarray[hash_vals[j0].second], myarray[hash_vals[j].second]))
      {
        gd = j;
        while (gd < j1)
        {
          int num = 0;
          for (size_t k = gd + 1; k < j1; k++)
          {
            if (pr(myarray[hash_vals[j0].second], myarray[hash_vals[k].second]))
            {
              num++;
              std::swap(hash_vals[k].second, hash_vals[gd].second);
              break;
            }
          }
          if (num > 0)
            gd++;
          else
          {
            j1 = gd;
            break;
          }
        }
        break;
      }
    }
    return j1;
  }


  template <typename ARRAY, typename Pr>
  void ArrayHashProxy<ARRAY, Pr>::cleanup(const DIntArray& data_temp, int cnt_new, ARRAY* deleted)
  {
    int i, gd0 = cnt_new + this->i0(), sz = this->_array.size();
    int ji, gd = gd0;
    ARRAY& marray = (ARRAY&)(this->_array);
    for (i = this->i0(); i < gd0; i++)
    {
      if (data_temp[i] > 0)
        continue;
      for (ji = gd; ji < sz; ji++)
      {
        if (data_temp[ji] > 0)
        {
          std::swap(marray[i], marray[ji]);
          gd = ji + 1;
          break;
        }
      }
    }

    if (deleted)
    {
      for (int i = cnt_new + this->i0(); i < sz; i++)
      {
        deleted->push_back(marray[i]);
      }
    }
    marray.resize(cnt_new + this->i0());
  }

  template <typename ARRAY, typename Pr>
  void ArrayHashProxy<ARRAY, Pr>::unique(ARRAY* deleted)
  {
    if (this->_array.empty())
      return;
    this->hash();
    size_t ji = this->i0(), sz = this->_array.size();
    int cnt_new = 0;
    DIntArray data_temp(sz, 0);
    while (ji < sz)
    {
      data_temp[(this->_hash_vals)[ji].second] = 1;
      ji = ArrayHashProxyC<ARRAY, Pr>::Clean(this->_array, this->_pr, (this->_hash_vals), ji, sz);
      ++cnt_new;
    }
    this->cleanup(data_temp, cnt_new, deleted);
    this->hash();
    (this->_clean) = 1;
  }

  template <typename ARRAY, typename Pr>
  void ArrayHashProxy<ARRAY, Pr>::erase_duplicate(ARRAY* deleted, int nErase)
  {
    if (this->_array.empty())
      return;
    this->hash();
    size_t ji = this->i0(), sz = this->_array.size();

    DIntArray data_temp(sz, 0);
    int cnt_new = 0, ecnt = 0;
    while (ji < sz)
    {
      size_t j0 = ji;
      ji = ArrayHashProxyC<ARRAY, Pr>::Clean(this->_array, this->_pr, (this->_hash_vals), ji, sz);
      int num = ji - j0;
      bool remainIt = false;
      if (nErase > 0)
      {
        if (num != nErase)
          remainIt = true;
      }
      else if (num == 1)
        remainIt = true;

      if (remainIt)
      {
        data_temp[(this->_hash_vals)[j0].second] = 1;
        cnt_new++;
      }
      else
        ecnt++;
    }
    if (ecnt < 1)
      return;
    this->cleanup(data_temp, cnt_new, deleted);
    (this->_clean) = 0;
    this->hash();
    this->_clean = 1;
  }

  template <typename ARRAY, typename Pr>
  void ArrayHashProxyC<ARRAY, Pr>::update()
  {
    if (this->_array.empty())
      return;
    this->hash();
    size_t ji = this->i0(), sz = this->_array.size();
    while (ji < sz)
    {
      ji = Clean(this->_array, this->_pr, (this->_hash_vals), ji, sz);
    }

    (this->_clean) = 1;
  }

  template <typename ARRAY, typename Pr>
  size_t ArrayHashProxyC<ARRAY, Pr>::hash_(const value_type& src)const
  {
    return this->_pr(src);
  }

  template <typename ARRAY, typename Pr>
  void ArrayHashProxyC<ARRAY, Pr>::hash()
  {
    if ((this->_hash_vals).size() != (this->i0() + this->_array.size()))
      (this->_clean) = 0;

    if (this->_array.empty() ||
      (this->_clean) == 1)
      return;
    int i, sz = this->_array.size();
    (this->_hash_vals).resize(sz);


    for (i = this->i0(); i < sz; i++)
    {
      (this->_hash_vals)[i].first = this->hash_(this->_array[i]);
      (this->_hash_vals)[i].second = i;
    }

    std::sort((this->_hash_vals).begin(), (this->_hash_vals).end());
  }

  inline
    bool SetByString(Json::Value& val, const std::string& str)
  {
    std::stringstream ss;
    ss << str;
    Json::Reader reader;
    return reader.parse(ss, val);
  }

  inline
    bool SetByFile(Json::Value& val, const std::string& fn)
  {
    std::ifstream fin(fn);
    if (!fin)
      return false;
    std::string fileStr;
    fileStr.insert(fileStr.begin(), std::istreambuf_iterator<char>(fin.rdbuf()), std::istreambuf_iterator<char>());
    return SetByString(val, fileStr);
  }

  inline
    bool SetToFile(const Json::Value& val, const std::string& fn)
  {
    std::ofstream fout(fn);
    if (!fout)
      return false;
    Json::FastWriter wr;
    auto fileStr = wr.write(val);
    auto len = strlen(fileStr.data());
    fout.write(fileStr.data(), len);
    fout.close();
    return true;
  }

  inline
    bool GetChild(const Json::Value& val, const std::string& str, Json::Value& rt)
  {
    try {
      auto pos = str.find_first_of('.');
      if (pos != std::string::npos)
      {
        Json::Value myChild;
        if (!GetChild(val, str.substr(0, pos), myChild))
          return false;
        return GetChild(myChild, str.substr(pos + 1), rt);
      }
      if (val.empty() || !val.isMember(str))
        return false;
      rt = val[str];
      return true;
    }
    catch (...) {
      return false;
    }
    return true;
  }

  inline
    bool GetValue(const Json::Value& val, const std::string& str, std::string& rt)
  {
    try {
      auto pos = str.find_first_of('.');
      if (pos != std::string::npos)
      {
        Json::Value nextValue;
        if (GetChild(val, str.substr(0, pos), nextValue))
        {
          return GetValue(nextValue, str.substr(pos + 1), rt);
        }
        else
          return false;
      }

      if (!val.isMember(str))
        return false;
      std::stringstream ss;
      ss << val[str].asString();
      rt = ss.str();
    }
    catch (...) {
      return false;
    }
    return true;
  }

  inline bool HasKey(const Json::Value& val, const std::string& str)
  {
      Json::Value child;
      return GetChild(val, str, child);
  }

  inline
    bool WriteValue(Json::Value& val, const std::string& str, const std::string& rt)
  {
    try {
      auto pos = str.find_first_of('.');
      if (pos != std::string::npos)
      {
        val[str.substr(0, pos)][str.substr(pos + 1)] = rt;
      }
      else
      {
        val[str] = rt;
      }
      return true;
    }
    catch (...) {
      return false;
    }
    return true;
  }

  template<typename ARRAY>
  bool GetData(const Json::Value& val, const std::string& str, ARRAY& data)
  {
    std::string strData, strOut;
    if (!GetValue(val, str, strData))
    {
      return false;
    }
    if (!Base64::Decode(strData, &strOut))
    {
      return false;
    }

    std::stringstream mIO;
    mIO.write(strOut.data(), strOut.size());

    int32_t num;

    mIO.read((char*)(&num), sizeof(int32_t));
    if (num < 1)
      return false;
    data.resize(num);
    mIO.read((char*)(data.data()), num * sizeof(typename ARRAY::value_type));
    return true;
  }


  template <typename T>
  bool tGetValue(const Json::Value &jsonData, const char *key, T &val)
  {
      try {
          Json::Value val0;
          if (!GetChild(jsonData, key, val0))
          {
              return false;
          }
          std::stringstream ss;
          if (val0.isUInt()) {
              ss << val0.asUInt();
          }
          else if (val0.isInt()) {
              ss << val0.asInt();
          }
          else if (val0.isBool()) {
              ss << val0.asBool();
          }
          else if (val0.isNumeric()) {
              ss << val0.asDouble();
          }
          else if (val0.isString()) {
              ss << val0.asString();
          }
          else {
              return false;
          }
          ss >> val;
      }
      catch (...) {
          return false;
      }
      return true;
  }

  inline
  bool GetData(const Json::Value& val, const std::string& str, 
    std::vector<double>& pointlist, std::vector<int>& facetmarkerlist, std::vector<int>& facetlist, std::vector<int>& tetrahedronlist, std::vector<int>& neighborlist)
  {

    std::string strData, strOut;
    if (!GetValue(val, str, strData))
    {
      return false;
    }

    if (!Base64::Decode(strData, &strOut))
    {
      return false;
    }

    std::stringstream mIO;
    mIO.write(strOut.data(), strOut.size());

    int numberoftetrahedra = 0, numberoffacets = 0, numberofpoints = 0;

    mIO.read((char*)(&numberoftetrahedra), sizeof(int));
    mIO.read((char*)(&numberoffacets), sizeof(int));
    mIO.read((char*)(&numberofpoints), sizeof(int));

    pointlist.resize(numberofpoints*3, 0);
    facetmarkerlist.resize(numberoffacets, 0);
    facetlist.resize(numberoffacets*3, 0);
    tetrahedronlist.resize(numberoftetrahedra*4, 0);
    neighborlist.resize(numberoftetrahedra * 4, 0);

    if (numberofpoints > 0) {
      mIO.read((char*)(pointlist.data()), sizeof(double)*numberofpoints * 3);
    }

    if (numberoffacets > 0)
    {
      mIO.read((char*)(facetmarkerlist.data()), sizeof(int)*numberoffacets);

      int numberofholes = 0, numberofpolygons = 1, numberofvertices = 3;
      for (int i = 0; i < numberoffacets; i++)
      {

        int f_pts[3];
        mIO.read((char*)(&numberofholes), sizeof(int));
        mIO.read((char*)(&numberofpolygons), sizeof(int));

        if (numberofholes > 0)
        {
          ;
        }

        if (numberofpolygons > 0)
        {
          for (int j = 0; j < numberofpolygons; j++) {
            mIO.read((char*)(&numberofvertices), sizeof(int));
            mIO.read((char*)(f_pts), sizeof(int)*numberofvertices);
            for (int k = 0; k < numberofvertices; k++)
              facetlist[i * 3 + k] = f_pts[k];
          }
        }

      }


    }

    if (numberoftetrahedra > 0)
    {
      mIO.read((char*)(tetrahedronlist.data()), sizeof(int)*numberoftetrahedra * 4);

      int numNeighbour = 0;
      mIO.read((char*)(&numNeighbour), sizeof(int));
      if (numNeighbour > 0) {
        mIO.read((char*)(neighborlist.data()), sizeof(int)*numberoftetrahedra * 4);
      }
    }

    return true;

  }

  inline
    bool GetDataByFile(const std::string& fn,
      std::vector<double>& pointlist, std::vector<int>& facetmarkerlist, std::vector<int>& facetlist, std::vector<int>& tetrahedronlist, std::vector<int>& neighborlist)
  {
    std::ifstream mIO(fn.c_str(), std::ios_base::in | std::ios_base::binary);

    int numberoftetrahedra = 0, numberoffacets = 0, numberofpoints = 0;

    mIO.read((char*)(&numberoftetrahedra), sizeof(int));
    mIO.read((char*)(&numberoffacets), sizeof(int));
    mIO.read((char*)(&numberofpoints), sizeof(int));

    pointlist.resize(numberofpoints * 3, 0);
    facetmarkerlist.resize(numberoffacets, 0);
    facetlist.resize(numberoffacets * 3, 0);
    tetrahedronlist.resize(numberoftetrahedra * 4, 0);
    neighborlist.resize(numberoftetrahedra * 4, 0);

    if (numberofpoints > 0) {
      mIO.read((char*)(pointlist.data()), sizeof(double)*numberofpoints * 3);
    }

    if (numberoffacets > 0)
    {
      mIO.read((char*)(facetmarkerlist.data()), sizeof(int)*numberoffacets);

      int numberofholes = 0, numberofpolygons = 1, numberofvertices = 3;
      for (int i = 0; i < numberoffacets; i++)
      {

        int f_pts[3];
        mIO.read((char*)(&numberofholes), sizeof(int));
        mIO.read((char*)(&numberofpolygons), sizeof(int));

        if (numberofholes > 0)
        {
          ;
        }

        if (numberofpolygons > 0)
        {
          for (int j = 0; j < numberofpolygons; j++) {
            mIO.read((char*)(&numberofvertices), sizeof(int));
            mIO.read((char*)(f_pts), sizeof(int)*numberofvertices);
            for (int k = 0; k < numberofvertices; k++)
              facetlist[i * 3 + k] = f_pts[k];
          }
        }

      }


    }

    if (numberoftetrahedra > 0)
    {
      mIO.read((char*)(tetrahedronlist.data()), sizeof(int)*numberoftetrahedra * 4);

      int numNeighbour = 0;
      mIO.read((char*)(&numNeighbour), sizeof(int));
      if (numNeighbour > 0) {
        mIO.read((char*)(neighborlist.data()), sizeof(int)*numberoftetrahedra * 4);
      }
    }

    mIO.close();
    return true;

  }

  inline
  bool WriteData(std::stringstream& mIO, std::string& strData)
  {
    std::string strIn;

    mIO.seekg(0, std::ios::end);
    int nLengthSize = mIO.tellg();
    mIO.seekg(0, std::ios::beg);

    strIn.resize(nLengthSize, 0);
    mIO.read((char*)strIn.data(), nLengthSize);


    if (!Base64::Encode(strIn, &strData))
    {
      return false;
    }

    return true;
  }


  template<typename ARRAY>
  bool WriteData(Json::Value& val, const std::string& str, const ARRAY& data)
  {
    std::stringstream mIO;
    std::string strData, strIn;

    int32_t num = data.size();

    if (num < 1)
      return false;

    mIO.write((char*)(&num), sizeof(int32_t));

    mIO.write((char*)(data.data()), num * sizeof(typename ARRAY::value_type));

    mIO.seekg(0, std::ios::end);
    int nLengthSize = mIO.tellg();
    mIO.seekg(0, std::ios::beg);

    strIn.resize(nLengthSize, 0);
    mIO.read((char*)strIn.data(), nLengthSize);


    if (!Base64::Encode(strIn, &strData))
    {
      return false;
    }

    if (!WriteValue(val, str, strData))
    {
      return false;
    }

    return true;
  }

  inline
  bool WriteData(Json::Value& val, const std::string& str, int numberoftetrahedra, int numberoffacets, int numberofpoints,
  double * pointlist, int  *facetmarkerlist, int *facetlist,  int *tetrahedronlist, int *neighborlist)
  {
    std::stringstream mIO;
    std::string strData;

    mIO.write((char*)(&numberoftetrahedra), sizeof(int));
    mIO.write((char*)(&numberoffacets), sizeof(int));
    mIO.write((char*)(&numberofpoints), sizeof(int));

    if (numberofpoints > 0) {
      mIO.write((char*)(pointlist), sizeof(double)*numberofpoints * 3);
    }

    if (numberoffacets > 0)
    {
      mIO.write((char*)(facetmarkerlist), sizeof(int)*numberoffacets);

      int numberofholes = 0, numberofpolygons = 1, numberofvertices = 3;
      for (int i = 0; i < numberoffacets; i++)
      {

        int f_pts[3] = { facetlist[i * 3], facetlist[i * 3 + 1], facetlist[i * 3 + 2] };
        mIO.write((char*)(&numberofholes), sizeof(int));
        mIO.write((char*)(&numberofpolygons), sizeof(int));

        if (numberofholes > 0)
        {
          ;
        }

        if (numberofpolygons > 0)
        {
          for (int j = 0; j < numberofpolygons; j++) {
            mIO.write((char*)(&numberofvertices), sizeof(int));
            mIO.write((char*)(f_pts), sizeof(int)*numberofvertices);
          }
        }
      }
    }

    if (numberoftetrahedra > 0)
    {
      mIO.write((char*)(tetrahedronlist), sizeof(int)*numberoftetrahedra * 4);

      int numNeighbour = neighborlist ? numberoftetrahedra : 0;
      mIO.write((char*)(&numNeighbour), sizeof(int));
      if (numNeighbour > 0) {
        mIO.write((char*)(neighborlist), sizeof(int)*numberoftetrahedra * 4);
      }
    }

    if (!Base64::Encode(mIO.str(), &strData))
    {
      return false;
    }

    if (!WriteValue(val, str, strData))
    {
      return false;
    }

    return true;

  }

  inline
    bool WriteData2File(Json::Value& val, const std::string& str, int numberoftetrahedra, int numberoffacets, int numberofpoints,
      double * pointlist, int  *facetmarkerlist, int *facetlist, int *tetrahedronlist, int *neighborlist)
  {
    std::string fn;
    GetValue(val, str, fn);
    std::fstream mIO(fn.c_str(), std::ios_base::out | std::ios_base::binary);

    std::string strData;

    mIO.write((char*)(&numberoftetrahedra), sizeof(int));
    mIO.write((char*)(&numberoffacets), sizeof(int));
    mIO.write((char*)(&numberofpoints), sizeof(int));

    if (numberofpoints > 0) {
      mIO.write((char*)(pointlist), sizeof(double)*numberofpoints * 3);
    }

    if (numberoffacets > 0)
    {
      mIO.write((char*)(facetmarkerlist), sizeof(int)*numberoffacets);

      int numberofholes = 0, numberofpolygons = 1, numberofvertices = 3;
      for (int i = 0; i < numberoffacets; i++)
      {

        int f_pts[3] = { facetlist[i * 3], facetlist[i * 3 + 1], facetlist[i * 3 + 2] };
        mIO.write((char*)(&numberofholes), sizeof(int));
        mIO.write((char*)(&numberofpolygons), sizeof(int));

        if (numberofholes > 0)
        {
          ;
        }

        if (numberofpolygons > 0)
        {
          for (int j = 0; j < numberofpolygons; j++) {
            mIO.write((char*)(&numberofvertices), sizeof(int));
            mIO.write((char*)(f_pts), sizeof(int)*numberofvertices);
          }
        }
      }
    }

    if (numberoftetrahedra > 0)
    {
      mIO.write((char*)(tetrahedronlist), sizeof(int)*numberoftetrahedra * 4);

      int numNeighbour = neighborlist ? numberoftetrahedra : 0;
      mIO.write((char*)(&numNeighbour), sizeof(int));
      if (numNeighbour > 0) {
        mIO.write((char*)(neighborlist), sizeof(int)*numberoftetrahedra * 4);
      }
    }

    mIO.close();

    return true;
  }


  void DebugBDF(std::ostream& out, size_t nnode, size_t nelem,
    std::function<bool(size_t i, int32_t &nid, double *pos)> nodeFun,
    std::function<bool(size_t i, int32_t &elemId, int32_t &group, int32_t &, int32_t *nodes)> elemFun);

  void DebugBDF(const std::string & fn, size_t nnode, size_t nelem,
    std::function<bool(size_t i, int32_t &, double *)> nodeFun,
    std::function<bool(size_t i, int32_t &, int32_t &, int32_t &, int32_t *)> elemFun);

};